# conda install -c conda-forge osmnx
# https://osmnx.readthedocs.io/en/stable/
# https://qiita.com/moshi/items/e383491664d028cd2166
#
# https://zenn.dev/homing/articles/f9a314841c737d
# https://qiita.com/naozo-se/items/bc168a2282907a6b29bf
# https://stackoverflow.com/questions/62285134/how-to-fill-water-bodies-with-osmnx-in-python
import osmnx as ox
import networkx as nx
import folium
from geopy.distance import distance
# 与えられた緯度経度
latitude = 34.69155489576708 # 例としてケーニッヒスベルクの緯度
longitude = 135.48959938567188 # 例としてケーニッヒスベルクの経度
# 2kmの範囲を定義
distance_km = 2
# 地図データを取得
G = ox.graph_from_point((latitude, longitude),
dist=distance_km*1000)
#network_type='drive')
# グラフを取得
nodes, edges = ox.graph_to_gdfs(G)
# 橋の情報を抽出('bridge'タグを持つエッジのみを抽出)
bridge_edges = edges[edges['bridge'] == 'yes']
#bridge_edges = edges
# グラフのノードとエッジを抽出
bridge_graph = nx.Graph()
# エッジ(橋)をグラフに追加
for idx, data in bridge_edges.iterrows():
u, v, key = idx
bridge_graph.add_edge(u, v)
# オイラー閉路の判定
eulerian = nx.is_eulerian(bridge_graph)
# Foliumで可視化
m = folium.Map(location=[latitude, longitude], zoom_start=14)
# グラフのエッジ(橋)を地図に追加
for u, v in bridge_graph.edges():
point_u = (G.nodes[u]['y'], G.nodes[u]['x'])
point_v = (G.nodes[v]['y'], G.nodes[v]['x'])
folium.PolyLine(locations=[point_u, point_v], color='blue').add_to(m)
# 中心地点にマーカーを追加
folium.Marker(location=[latitude, longitude], popup="中心地点").add_to(m)
# オイラー閉路の有無を表示
if eulerian:
folium.Marker(location=[latitude, longitude], popup="オイラー閉路が存在します。", icon=folium.Icon(color='green')).add_to(m)
else:
folium.Marker(location=[latitude, longitude], popup="オイラー閉路は存在しません。", icon=folium.Icon(color='red')).add_to(m)
m
Make this Notebook Trusted to load map: File -> Trust Notebook
import osmnx as ox
import networkx as nx
import folium
# 与えられた緯度経度
latitude = 54.7088 # 例としてケーニッヒスベルクの緯度
longitude = 20.5108 # 例としてケーニッヒスベルクの経度
# 2kmの範囲を定義
distance_km = 2
# 地図データを取得
G = ox.graph_from_point((latitude, longitude), dist=distance_km*1000, network_type='all')
# 無向グラフに変換
G = G.to_undirected()
# グラフの連結成分を取得
connected_components = list(nx.connected_components(G))
# 橋の情報を抽出('bridge'タグを持つエッジのみを抽出)
bridge_edges = [edge for edge in G.edges(data=True) if 'bridge' in edge[2] and edge[2]['bridge'] == 'yes']
# 各連結成分を一つのノード(陸地)と見なす
land_nodes = {i: component for i, component in enumerate(connected_components)}
# 橋のノードを含む陸地ノードのIDを見つける関数
def find_land_node(node):
for land_id, nodes in land_nodes.items():
if node in nodes:
return land_id
return None
# 橋で陸地ノードを接続するための新しいグラフを作成
final_graph = nx.Graph()
# 陸地ノードを追加
for land_id, nodes in land_nodes.items():
final_graph.add_node(land_id, nodes=nodes)
# 橋エッジを追加(橋が異なる陸地を接続している場合)
for u, v, data in bridge_edges:
land_u = find_land_node(u)
land_v = find_land_node(v)
if land_u is not None and land_v is not None and land_u != land_v:
final_graph.add_edge(land_u, land_v, **data)
# オイラー閉路の判定
eulerian = nx.is_eulerian(final_graph)
# Foliumで可視化
m = folium.Map(location=[latitude, longitude], zoom_start=14)
# 橋のエッジを地図に追加
for u, v, data in final_graph.edges(data=True):
land_u_nodes = land_nodes[u]
land_v_nodes = land_nodes[v]
point_u = next((G.nodes[node]['y'], G.nodes[node]['x']) for node in land_u_nodes)
point_v = next((G.nodes[node]['y'], G.nodes[node]['x']) for node in land_v_nodes)
folium.PolyLine(locations=[point_u, point_v], color='blue').add_to(m)
# 中心地点にマーカーを追加
folium.Marker(location=[latitude, longitude], popup="中心地点").add_to(m)
# オイラー閉路の有無を表示
if eulerian:
folium.Marker(location=[latitude, longitude], popup="オイラー閉路が存在します。", icon=folium.Icon(color='green')).add_to(m)
else:
folium.Marker(location=[latitude, longitude], popup="オイラー閉路は存在しません。", icon=folium.Icon(color='red')).add_to(m)
m
Make this Notebook Trusted to load map: File -> Trust Notebook
print(land_graph)
print(final_graph)
MultiGraph with 5819 nodes and 8690 edges
Graph with 1 nodes and 0 edges
import folium
import osmnx as ox
import geopandas as gpd
# 対象地域の道路情報取得 (愛知県名古屋市中村区)
#query = "Nakamuraku,Nagoya,Aichi,Japan"
#G = ox.graph_from_place(query, network_type="drive")
(lat,lon)=(34.69199639822432, 135.49084584150586)
# Bounding box
bN,bS,bE,bW=lon+1,lon-1,lon-1,lon+1
#graph=ox.graph.graph_from_point((lat, lon),
# dist=1000,
# network_type='drive')
graph=ox.geometries.geometries_from_bbox(bN, bS, bE, bW,
tags={"natural": "water"})
# 道路グラフネットワーク可視化
#fmap=ox.plot_graph_folium(graph)
fmap=ox.plot_graph_folium(graph)
folium.LayerControl().add_to(fmap)
fmap
#fmap.save(outfile="road_network.html")
#opts = {"node_size": 5, "bgcolor": "white", "node_color": "blue", "edge_color": "blue"}
#ox.plot_graph(G, show=False, save=True, filepath="road_network.png", **opts)
/var/folders/mn/nmxx09wj21j7f81mwghvt1sw0000gn/T/ipykernel_40987/1645265506.py:18: FutureWarning: The `geometries` module and `geometries_from_X` functions have been renamed the `features` module and `features_from_X` functions. Use these instead. The `geometries` module and function names are deprecated and will be removed in the v2.0.0 release. See the OSMnx v2 migration guide: https://github.com/gboeing/osmnx/issues/1123
graph=ox.geometries.geometries_from_bbox(bN, bS, bE, bW,
/opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/geometries.py:48: FutureWarning: The `north`, `south`, `east`, and `west` parameters are deprecated and will be removed in the v2.0.0 release. Use the `bbox` parameter instead. See the OSMnx v2 migration guide: https://github.com/gboeing/osmnx/issues/1123
return features.features_from_bbox(north, south, east, west, tags=tags)
---------------------------------------------------------------------------
IndexError Traceback (most recent call last)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/geopandas/array.py:950, in GeometryArray.estimate_utm_crs(self, datum_name)
949 try:
--> 950 return CRS.from_epsg(utm_crs_list[0].code)
951 except IndexError:
IndexError: list index out of range
During handling of the above exception, another exception occurred:
RuntimeError Traceback (most recent call last)
Cell In[17], line 18
13 bN,bS,bE,bW=lon+1,lon-1,lon-1,lon+1
15 #graph=ox.graph.graph_from_point((lat, lon),
16 # dist=1000,
17 # network_type='drive')
---> 18 graph=ox.geometries.geometries_from_bbox(bN, bS, bE, bW,
19 tags={"natural": "water"})
21 # 道路グラフネットワーク可視化
22 #fmap=ox.plot_graph_folium(graph)
23 fmap=ox.plot_graph_folium(graph)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/geometries.py:48, in geometries_from_bbox(north, south, east, west, tags)
22 """
23 Do not use: deprecated.
24
(...)
45 gdf : geopandas.GeoDataFrame
46 """
47 warn(DEP_MSG, FutureWarning, stacklevel=2)
---> 48 return features.features_from_bbox(north, south, east, west, tags=tags)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:134, in features_from_bbox(north, south, east, west, bbox, tags)
131 polygon = utils_geo.bbox_to_poly(bbox=bbox)
133 # create GeoDataFrame of features within this polygon
--> 134 return features_from_polygon(polygon, tags)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:349, in features_from_polygon(polygon, tags)
346 response_jsons = _overpass._download_overpass_features(polygon, tags)
348 # create GeoDataFrame from the downloaded data
--> 349 return _create_gdf(response_jsons, polygon, tags)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:440, in _create_gdf(response_jsons, polygon, tags)
437 relation_types = {"boundary", "multipolygon"}
439 # extract geometries from the downloaded osm data
--> 440 for response_json in response_jsons:
441 response_count += 1
443 # Parses the JSON of OSM nodes, ways and (multipolygon) relations
444 # to dictionaries of coordinates, Shapely Points, LineStrings,
445 # Polygons and MultiPolygons
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/_overpass.py:389, in _download_overpass_features(polygon, tags)
373 """
374 Retrieve OSM features within boundary from the Overpass API.
375
(...)
386 a generator of JSON responses from the Overpass server
387 """
388 # subdivide query polygon to get list of sub-divided polygon coord strings
--> 389 polygon_coord_strs = _make_overpass_polygon_coord_strs(polygon)
390 utils.log(f"Requesting data from API in {len(polygon_coord_strs)} request(s)")
392 # pass exterior coordinates of each polygon in list to API, one at a time
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/_overpass.py:253, in _make_overpass_polygon_coord_strs(polygon)
234 """
235 Subdivide query polygon and return list of coordinate strings.
236
(...)
249 list of strings of exterior coordinates of polygon(s)
250 """
251 # first subdivide the polygon if its area exceeds max size
252 # this results in a multipolygon of 1+ constituent polygons
--> 253 poly_proj, crs_proj = projection.project_geometry(polygon)
254 multi_poly_proj = utils_geo._consolidate_subdivide_geometry(poly_proj)
255 multi_poly, _ = projection.project_geometry(multi_poly_proj, crs=crs_proj, to_latlong=True)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/projection.py:60, in project_geometry(geometry, crs, to_crs, to_latlong)
57 crs = settings.default_crs
59 gdf = gpd.GeoDataFrame(geometry=[geometry], crs=crs)
---> 60 gdf_proj = project_gdf(gdf, to_crs=to_crs, to_latlong=to_latlong)
61 geometry_proj = gdf_proj["geometry"].iloc[0]
62 return geometry_proj, gdf_proj.crs
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/projection.py:99, in project_gdf(gdf, to_crs, to_latlong)
97 # else if to_crs is None, project gdf to an appropriate UTM zone
98 elif to_crs is None:
---> 99 to_crs = gdf.estimate_utm_crs()
101 # project the gdf
102 gdf_proj = gdf.to_crs(to_crs)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/geopandas/geodataframe.py:1451, in GeoDataFrame.estimate_utm_crs(self, datum_name)
1415 def estimate_utm_crs(self, datum_name="WGS 84"):
1416 """Returns the estimated UTM CRS based on the bounds of the dataset.
1417
1418 .. versionadded:: 0.9
(...)
1449 - Prime Meridian: Greenwich
1450 """
-> 1451 return self.geometry.estimate_utm_crs(datum_name=datum_name)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/geopandas/geoseries.py:1214, in GeoSeries.estimate_utm_crs(self, datum_name)
1178 def estimate_utm_crs(self, datum_name: str = "WGS 84") -> CRS:
1179 """Returns the estimated UTM CRS based on the bounds of the dataset.
1180
1181 .. versionadded:: 0.9
(...)
1212 - Prime Meridian: Greenwich
1213 """
-> 1214 return self.values.estimate_utm_crs(datum_name)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/geopandas/array.py:952, in GeometryArray.estimate_utm_crs(self, datum_name)
950 return CRS.from_epsg(utm_crs_list[0].code)
951 except IndexError:
--> 952 raise RuntimeError("Unable to determine UTM CRS")
RuntimeError: Unable to determine UTM CRS
import osmnx as ox
ox.settings.log_console = True
# add more items here to get all the landforms you want
places = ['Manhattan, NY, USA', 'Brooklyn, NY, USA', 'Queens, NY, USA', 'Bronx, NY, USA']
land = ox.geocode_to_gdf(places)
# get the water bodies
left, bottom, right, top = land.total_bounds
bbox = top, bottom, right, left
poly = ox.utils_geo.bbox_to_poly(*bbox)
water = ox.features_from_polygon(poly, tags={'natural': 'water'})
# constrain the plotting window as desired
c = land.unary_union.centroid
bbox = ox.utils_geo.bbox_from_point((c.y, c.x), dist=12000)
water_color = 'blue'
land_color = '#aaaaaa'
fig, ax = ox.plot_footprints(water, bbox=bbox,
color=water_color, bgcolor=water_color,
show=False, close=False)
ax = land.plot(ax=ax, zorder=0, fc=land_color)
/var/folders/mn/nmxx09wj21j7f81mwghvt1sw0000gn/T/ipykernel_40987/273932149.py:11: FutureWarning: The `north`, `south`, `east`, and `west` parameters are deprecated and will be removed in the v2.0.0 release. Use the `bbox` parameter instead. See the OSMnx v2 migration guide: https://github.com/gboeing/osmnx/issues/1123
poly = ox.utils_geo.bbox_to_poly(*bbox)
import osmnx as ox
import folium
# 場所を設定(例:東京都渋谷区)
place_name = "Shibuya, Tokyo, Japan"
(lat,lon)=(34.69199639822432, 135.49084584150586)
# 陸地と道路のグラフを取得
#G = ox.graph_from_place(place_name, network_type='drive')
G = ox.graph_from_point((lat,lon), network_type='drive')
nodes, edges = ox.graph_to_gdfs(G)
# 中心点を取得
lat = nodes['y'].mean()
lon = nodes['x'].mean()
# foliumマップを作成
m = folium.Map(location=[lat, lon], zoom_start=14)
# 陸地を描画
folium.GeoJson(data=edges).add_to(m)
# マップを保存または表示
m.save('map.html')
# マップを表示(Jupyter Notebookの場合)
m
Make this Notebook Trusted to load map: File -> Trust Notebook
import osmnx as ox
import folium
import geopandas as gpd
# 場所を設定(例:東京都渋谷区)
place_name = "Shibuya, Tokyo, Japan"
# 陸地のポリゴンを取得
land_polygons = ox.geometries.geometries_from_place(place_name, tags={'natural': 'land'})
# 水域のポリゴンを取得
water_polygons = ox.geometries.geometries_from_place(place_name, tags={'natural': ['water', 'bay', 'strait']})
# 道路のグラフを取得
G = ox.graph_from_place(place_name, network_type='drive')
nodes, edges = ox.graph_to_gdfs(G)
# 中心点を取得
lat = nodes['y'].mean()
lon = nodes['x'].mean()
# foliumマップを作成
m = folium.Map(location=[lat, lon], zoom_start=14)
# 陸地を描画
for _, polygon in land_polygons.iterrows():
folium.GeoJson(polygon['geometry'], style_function=lambda x: {'fillColor': 'green', 'color': 'green'}).add_to(m)
# 水域を描画
for _, polygon in water_polygons.iterrows():
folium.GeoJson(polygon['geometry'], style_function=lambda x: {'fillColor': 'blue', 'color': 'blue'}).add_to(m)
# 道路を描画
folium.GeoJson(data=edges).add_to(m)
# マップを保存または表示
m.save('map.html')
# マップを表示(Jupyter Notebookの場合)
m
/var/folders/mn/nmxx09wj21j7f81mwghvt1sw0000gn/T/ipykernel_40987/388463500.py:9: FutureWarning: The `geometries` module and `geometries_from_X` functions have been renamed the `features` module and `features_from_X` functions. Use these instead. The `geometries` module and function names are deprecated and will be removed in the v2.0.0 release. See the OSMnx v2 migration guide: https://github.com/gboeing/osmnx/issues/1123
land_polygons = ox.geometries.geometries_from_place(place_name, tags={'natural': 'land'})
---------------------------------------------------------------------------
InsufficientResponseError Traceback (most recent call last)
Cell In[22], line 9
6 place_name = "Shibuya, Tokyo, Japan"
8 # 陸地のポリゴンを取得
----> 9 land_polygons = ox.geometries.geometries_from_place(place_name, tags={'natural': 'land'})
11 # 水域のポリゴンを取得
12 water_polygons = ox.geometries.geometries_from_place(place_name, tags={'natural': ['water', 'bay', 'strait']})
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/geometries.py:128, in geometries_from_place(query, tags, which_result, buffer_dist)
104 """
105 Do not use: deprecated.
106
(...)
125 gdf : geopandas.GeoDataFrame
126 """
127 warn(DEP_MSG, FutureWarning, stacklevel=2)
--> 128 return features.features_from_place(query, tags, which_result, buffer_dist)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:298, in features_from_place(query, tags, which_result, buffer_dist)
295 utils.log("Constructed place geometry polygon(s) to query API")
297 # create GeoDataFrame using this polygon(s) geometry
--> 298 return features_from_polygon(polygon, tags)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:349, in features_from_polygon(polygon, tags)
346 response_jsons = _overpass._download_overpass_features(polygon, tags)
348 # create GeoDataFrame from the downloaded data
--> 349 return _create_gdf(response_jsons, polygon, tags)
File /opt/anaconda3/envs/jupyter-book-py311/lib/python3.11/site-packages/osmnx/features.py:489, in _create_gdf(response_jsons, polygon, tags)
487 if len(geometries) == 0: # pragma: no cover
488 msg = "No data elements in server response. Check log and query location/tags."
--> 489 raise InsufficientResponseError(msg)
491 # remove untagged elements from the final dict of geometries
492 utils.log(f"{len(geometries)} geometries created in the dict")
InsufficientResponseError: No data elements in server response. Check log and query location/tags.
# 試しコード
import osmnx as ox
import networkx as nx
import folium
# 地図データを取得
# c.f. https://osmnx.readthedocs.io/en/stable/user-reference.html#module-osmnx.graph
graph=ox.graph.graph_from_bbox(
# bbox : (north, south, east, west)
bbox=(35.05,34.90,
135.75,135.80),
# "all", "all_public", "bike", "drive", "drive_service", "walk"
network_type='bike',
simplify=True, # True
retain_all=False,
truncate_by_edge=False,
custom_filter=None)
# 無向グラフにしておく
graph=graph.to_undirected()
# グラフを取得
nodes, edges = ox.graph_to_gdfs(graph)
# グラフのノードとエッジを追加
nxgraph = nx.Graph()
for idx, data in edges.iterrows():
u, v, key = idx
nxgraph.add_edge(u, v)
# Foliumで可視化
lat, lon = [(35.05+34.90)/2,
(135.80+135.75)/2]
fmap = folium.Map(location=[latitude, longitude], zoom_start=13)
# グラフを地図に追加
for u, v in nxgraph.edges():
point_u = (graph.nodes[u]['y'],graph.nodes[u]['x'])
point_v = (graph.nodes[v]['y'],graph.nodes[v]['x'])
folium.PolyLine(locations=[point_u, point_v], color='blue').add_to(fmap)
fmap
Make this Notebook Trusted to load map: File -> Trust Notebook